FIELDimageR:
Performing image analysis in plant breeding programs
Introduction
Image phenotyping platforms
Image phenotyping platforms are a fast and non-invasive tool for phenotyping plant under lab and field conditions. The most common are sensors for imaging. The reflected light captured in these images can be used to draw inference about many traits, including:
- Geometric traits (i.e. plant height, leaf area index, lodging, crop canopy cover)
- Canopy spectral texture (spectral features)
- Physiological traits (i.e., chlorophyll, biomass, pigment content, photosynthesis)
- Abiotic/biotic stress indicators (i.e., stomatal conductance, canopy temperature difference, leaf water potential, senescence index)
- Nutrients (nitrogen concentration, protein content)
- Yield
Image segmentation
Image segmentation is the process of partitioning a digital image into multiple segments (sets of pixels or image objects).
Thresholding
Thresholding is the simplest method of image segmentation. This method is based on a clip-level (or a threshold value) to turn a gray-scale image into a binary image (Mask).
FIELDimageR
FIELDimageR is a R package to analyze images from research fields, and allows to:
- Crop the image
- Remove the background
- Build vegetation indices
- Rotate the image
- Build the plot shapefile
- Extract information for each plot
1) Example 01
Counting green seeds in Soybean
Alexandre Abdão(Embrapa Maize & Sorghum)
4 pictures: soybean samples with different green seeds percentage
Number of green seeds (0, 1, 6, and 11)
Donwload:
soybean.zip
Required packages
# Installing
install.packages("sp")
install.packages("raster")
install.packages("rgdal")
install.packages("ggplot2")
install.packages("agricolae")
install.packages("reshape2")
install.packages("devtools")
install.packages("parallel")
install.packages("foreach")
install.packages("doParallel")
install.packages("lme4")
devtools::install_github("filipematias23/FIELDimageR")
# Necessary packages
library(FIELDimageR)
library(raster)
library(ggplot2)
library(agricolae)
library(reshape2)
library(lme4)- Uploading one image as example and decreasing the resolution
# Select one index to identify leaves and remove the background
EX1.I1<- fieldIndex(mosaic = EX1,index = c("SI","BGI","BI"))# Removing the background
EX1.R1<- fieldMask(mosaic = EX1, index = "BGI",
cropValue = 0.7,
cropAbove = T)# Counting the total number of seeds
EX.P.Total<-fieldCount(mosaic = EX1.R1$mask,
fieldShape = EX.shapeFile$fieldShape,
minSize = 0.1,
cex = 1.5,
na.rm = T)# Select one index to identify green seeds
EX1.I2<- fieldIndex(mosaic = EX1.R1$newMosaic,index = c("SI","BGI","BI"))# Selecting green seeds
EX1.R2<- fieldMask(mosaic = EX1.R1$newMosaic,
index = "BI",
#myIndex = "Blue",
cropValue = 130,
cropAbove = T)# Counting the number of green seeds
EX.Green<-fieldCount(mosaic = EX1.R2$mask,
fieldShape = EX.shapeFile$fieldShape,
minSize = 0.07,
cex = 1.5,
na.rm = T)# Joying information
data.frame(Total=EX.P.Total$fieldCount,
Green=EX.Green$fieldCount,
Percentage=round(EX.Green$fieldCount/EX.P.Total$fieldCount,2))
################
### Parallel ###
################
# Required packages
library(parallel)
library(foreach)
library(doParallel)
# Images names (folder directory: "./soybean/")
pics<-list.files("./soybean/")
# Number of cores
n.core<-2
# Starting parallel
cl <- makeCluster(n.core, output = "")
registerDoParallel(cl)
EX.Table.Parallel <- foreach(i = 1:length(pics), .packages = c("raster","FIELDimageR"),
.combine = rbind) %dopar% {
EX1<-stack(paste("./soybean/",pics[i],sep = ""))
EX1<-aggregate(EX1, fact= 4)
EX.shapeFile<-fieldPolygon(EX1,extent = T,
plot = F)
EX1.R1<- fieldMask(mosaic = EX1, index = "BGI",
cropValue = 0.7,
cropAbove = T,
plot = F)
EX.P.Total<-fieldCount(mosaic = EX1.R1$mask,
fieldShape = EX.shapeFile$fieldShape,
minSize = 0.1,
cex = 1.5,
na.rm = T)
EX1.R2<- fieldMask(mosaic = EX1.R1$newMosaic,
index = "BI",
cropValue = 130,
cropAbove = T,
plot = F)
EX.Green<-fieldCount(mosaic = EX1.R2$mask,
fieldShape = EX.shapeFile$fieldShape,
minSize = 0.07,
cex = 1.5,
na.rm = T)
data.frame(Total=EX.P.Total$fieldCount,
Green=EX.Green$fieldCount,
Percentage=round(EX.Green$fieldCount/EX.P.Total$fieldCount,2))
}
rownames(EX.Table.Parallel)<-pics
EX.Table.Parallel
2) Example 02
Remote sensing (Potato Breeding)
Jeffrey Endelman and Filipe Matias (UW-Madison)
- Donwload:
odm_orthophoto.tif
- Evaluating each trial per time:
- Dataset - Field notes - phenotype (Download: ‘EX1_Data.xlsx’)
# Dataset - Field notes - phenotype (.exel or .txt):
Data1<-read_excel('EX1_Data.xlsx',1)
# Field map ID identification (ID="Plot"):
Map1<-fieldMap(fieldPlot = Data1$Plot,fieldColumn = Data1$Column, fieldRow = Data1$Row,decreasing = T)
Map1# Building the plot shapefile (ncols = 14 and nrows = 10)
x11()
EX1.Shape.1<-fieldShape(mosaic = EX1.Crop.1, ncols = 14, nrows = 10, fieldData = Data1, ID = "Plot", fieldMap = Map1)# Trial: 02
EX1.Crop.2 <- fieldCrop(mosaic = EX1.RemSoil$newMosaic)
Data2<-read_excel('EX1_Data.xlsx',2)
Map2<-fieldMap(fieldPlot = Data2$Plot,fieldColumn = Data2$Column, fieldRow = Data2$Row,decreasing = T)
EX1.Shape.2<-fieldShape(mosaic = EX1.Crop.2, ncols = 14, nrows = 9, fieldData = Data2, ID = "Plot", fieldMap = Map2)# Trial: 03
EX1.Crop.3 <- fieldCrop(mosaic = EX1.RemSoil$newMosaic)
Data3<-read_excel('EX1_Data.xlsx',3)
Map3<-fieldMap(fieldPlot = Data3$Plot,fieldColumn = Data3$Column, fieldRow = Data3$Row,decreasing = T)
EX1.Shape.3<-fieldShape(mosaic = EX1.Crop.3, ncols = 14, nrows = 13, fieldData = Data3, ID = "Plot", fieldMap = Map3)# Trial: 04
EX1.Crop.4 <- fieldCrop(mosaic = EX1.RemSoil$newMosaic)
Data4<-read_excel('EX1_Data.xlsx',4)
Map4<-fieldMap(fieldPlot = Data4$Plot,fieldColumn = Data4$Column, fieldRow = Data4$Row,decreasing = T)
EX1.Shape.4<-fieldShape(mosaic = EX1.Crop.4, ncols = 14, nrows = 10, fieldData = Data4, ID = "Plot", fieldMap = Map4)# Combining shapefiles
EX1.Shape<-rbind(EX1.Shape.1$fieldShape,EX1.Shape.2$fieldShape,EX1.Shape.3$fieldShape,EX1.Shape.4$fieldShape)
plot(EX1.Shape)# Building indices ("NGRDI","BGI","GLI","SCI")
EX1.Indices<- fieldIndex(mosaic = EX1.RemSoil$newMosaic,
index = c("NGRDI","BGI"),
myIndex = c("(Red-Blue)/Green"))# Extracting Data
EX1.Info<- fieldInfo(mosaic = EX1.Indices[[c("NGRDI","BGI","myIndex")]],
fieldShape = EX1.Shape,
buffer = -0.05,
n.core = 3)
NewData<-EX1.Info$fieldShape@data
NewData# Making map plots
fieldPlot(fieldShape=EX1.Info$fieldShape,
fieldAttribute="NGRDI",
mosaic=EX1.Indices,
color=c("red","blue"),
alpha = 0.4,
round = 2)
- Estimation of plant height using a high throughput phenotyping platform. (Donwload:
dsm.tif)
# Extracting the estimate plant height average (EPH):
EPH <- fieldInfo(DSM.S$newMosaic, fieldShape = EX1.Info$fieldShape, fun = "mean",n.core = 3)
EPH$plotValue# Correlation
NewData2<-EPH$fieldShape@data[,c("layer","NGRDI","BGI","myIndex")]
cor1<-correlation(NewData2)
cor1$correlation
cor1$pvalue# Regression
NewData3<-EPH$fieldShape@data[,c("Trial","Name","Block","Row","Column","layer","NGRDI")]
NewData3$Trial<-as.factor(NewData3$Trial)
NewData3$Name<-as.factor(NewData3$Name)
NewData3$Block<-as.factor(NewData3$Block)
NewData3$Row<-as.factor(NewData3$Row)
NewData3$Column<-as.factor(NewData3$Column)
NewData3$layer<-as.numeric(as.character(NewData3$layer))
NewData3$NGRDI<-as.numeric(as.character(NewData3$NGRDI))
ggplot(NewData3,aes(x=layer,y=NGRDI, col=Trial))+
facet_wrap(~Trial, scales = "fixed",ncol = 1)+
geom_point() +
geom_smooth(method=lm)+
theme_bw()fieldPlot(fieldShape=EPH$fieldShape,
fieldAttribute="layer",
mosaic=EX1.Indices,
color=c("red","blue"),
alpha = 0.4,
round = 2)
- Heritability
# NGRDI
str(NewData3)
mod.T1<-lmer(NGRDI~Row+Column+(1|Name),NewData3[NewData3$Trial=="T1",])
H2.T1<-as.data.frame(VarCorr(mod.T1))$vcov[1]/sum(as.data.frame(VarCorr(mod.T1))$vcov)
mod.T2<-lmer(NGRDI~Row+Column+(1|Name),NewData3[NewData3$Trial=="T2",])
H2.T2<-as.data.frame(VarCorr(mod.T2))$vcov[1]/sum(as.data.frame(VarCorr(mod.T2))$vcov)
mod.T3<-lmer(NGRDI~Row+Column+(1|Name),NewData3[NewData3$Trial=="T3",])
H2.T3<-as.data.frame(VarCorr(mod.T3))$vcov[1]/sum(as.data.frame(VarCorr(mod.T3))$vcov)
mod.T4<-lmer(NGRDI~Row+Column+(1|Name),NewData3[NewData3$Trial=="T4",])
H2.T4<-as.data.frame(VarCorr(mod.T4))$vcov[1]/sum(as.data.frame(VarCorr(mod.T4))$vcov)
NewData4<-data.frame(Trail=factor(c("T1","T2","T3","T4")),
H2=as.numeric(c(H2.T1,H2.T2,H2.T3,H2.T4)))
ggplot(NewData4,aes(x=Trail,y=H2,fill=Trail))+
geom_bar(stat="identity")+
geom_text(aes(label=round(H2,2)), vjust=1.6, color="white", size=6)+
theme_bw()# Estimated plant height
mod.T1<-lmer(layer~Row+Column+(1|Name),NewData3[NewData3$Trial=="T1",])
H2.T1<-as.data.frame(VarCorr(mod.T1))$vcov[1]/sum(as.data.frame(VarCorr(mod.T1))$vcov)
mod.T2<-lmer(layer~Row+Column+(1|Name),NewData3[NewData3$Trial=="T2",])
H2.T2<-as.data.frame(VarCorr(mod.T2))$vcov[1]/sum(as.data.frame(VarCorr(mod.T2))$vcov)
mod.T3<-lmer(layer~Row+Column+(1|Name),NewData3[NewData3$Trial=="T3",])
H2.T3<-as.data.frame(VarCorr(mod.T3))$vcov[1]/sum(as.data.frame(VarCorr(mod.T3))$vcov)
mod.T4<-lmer(layer~Row+Column+(1|Name),NewData3[NewData3$Trial=="T4",])
H2.T4<-as.data.frame(VarCorr(mod.T4))$vcov[1]/sum(as.data.frame(VarCorr(mod.T4))$vcov)
NewData4<-data.frame(Trail=factor(c("T1","T2","T3","T4")),
EPH=as.numeric(c(H2.T1,H2.T2,H2.T3,H2.T4)))
ggplot(NewData4,aes(x=Trail,y=EPH,fill=Trail))+
geom_bar(stat="identity")+
geom_text(aes(label=round(EPH,2)), vjust=1.6, color="white", size=6)+
theme_bw()
Reference
Matias, F.I. (2019). FIELDimageR Pipeline (tutorial). https://github.com/filipematias23/FIELDimageR
Matias, F.I.; Caraza-Harter, M.; Endelman, J.B. (2020). FIELDimageR: A R Package to Analyze Orthomosaic Images from Agricultural Field Trials. The Plant Phenome Journal. DOI:10.1002/ppj2.20005